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Numerical simulations of turbulent channel flows, with or without additives, are limited in the 
extent of the Reynolds number He and Deborah number De. The comparison of such simulations to 
theories of drag reduction, which are usually derived for asymptotically high 7& and calls for 
some care. In this paper we present a study of drag reduction by rodlike polymers in a turbulent 
channel flow using direct numerical simulation and illustrate how these numerical results should be 
related to the recently developed theory. 

PACS numbers: xxxx 



INTRODUCTION 

Drag reduction in wall-bounded turbulent flows can 
be achieved with the addition of either flexible or rodlike 
polymers [l], Q . The theory of drag reduction by either 
type of polymers is by now well established [3]. The 
theory for flexible polymers is supported by both exper- 
iments and numerical simulations while that for rodlike 
polymers has been compared mainly with laboratory ex- 
periments. More scant are numerical simulations of drag 
reduction by rodlike polymers, and some results became 
available only quite recently [4]. Comparisons of the- 
ory and simulations in this case must be done with care, 
since the available theory assumes high enough Reynolds 
7fe and Deborah numbers in contrast to the situation 
in simulations where these crucial characteristic numbers 
are relatively low. To clarify the relation between the the- 
ory and the available numerical simulations, we present 
in this note a comparison between numerical simulations 
and theoretical predictions adapted to the limited values 
of 7fe and Ve . 

The equations of motion used for the numerical simu- 
lations are 
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supplemented by the incompressibility constraint 
dUi/dxi = where U is the velocity field, p is the 
pressure, the units are chosen such that the fluid density 
is unity, v$ is the kinematic viscosity of the neat fluid, 
and <7ij is the additional stress tensor due to the rodlike 
polymers. 

The rodlike polymers are represented by rigid and neu- 
trally buoyant elongated particles. The particles are as- 
sumed to be massless and have no inertia. The orienta- 
tion of each polymer is given by a unit vector n. In turbu- 
lent flows with strong shear such that thermal Brownian 



rotations can be neglected, the evolution equation for the 
second moment of the conformation tensor IZij = (niUj) 
is given by 

' ' ^fc 0^ = SikTZ-kj + SjkTlki — 2Ski7Zijki (2) 



dt 



dx k 



where TZijki = (niUjUkni) and Sij are the velocity gra- 
dient Sij = dui/dxj. Moreover, in this non-Brownian 
limit, the tensor cr^- in Eq. (pQ) is given by [1,0], 



<7ij = GvpSkiTZijki , 



(3) 



where v v is the polymeric contribution to the viscosity 
at vanishingly small and time-independent shear, and is 
proportional to the product of v§(\) where qb is the volume 
fraction of the polymers. In [5[ it was argued that these 
equations can be closed by the simple closure IZijki = 



SIMULATIONS 

The equations of motion are numerically integrated for 
a channel flow . The dimensions of the integration do- 
main are 2ttL x2Lx 1.2ttL in the x (streamwise), y (wall- 
normal) and z (spanwise) directions respectively, with L 
being the channel half-width. The numerical formula- 
tion is a standard pseudospectral method with Fourier 
expansion in directions parallel to the wall and Cheby- 
shev in the direction normal to the wall. The grid used is 
128 x 193 x 64. The direct numerical simulations (DNS) 
were performed at a nominal TZe = UqL/vq of 10000 for 
both the Newtonian flow (with r] p = 6v p /vo = 0) and 
for turbulent flow with rodlike polymers with r] p = 25, 
where Uq is the mean velocity at the center of the chan- 
nel. Both flows were forced on average with the same 
pressure drop p f = —dp/dx, so the resulting Reynolds 
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FIG. 1: Mean streamwise velocity profile V + (y + ) for both 
the Newtonian flow (solid line) and the rodlike polymer laden 
flow (circles). 
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FIG. 2: The Reynolds stress W(y) for the Newtonian flow 
(dashed line) and the flow with rodlike polymers (squares). 



COMPARISON OF THEORY WITH 
SIMULATIONS 



number based on the friction velocity is the same and 
equal to 7fe T = ^p'LL/vq = 300. In channel geome- 
try the only non-vanishing mean velocity component is 
V(y) = (U x ). Accordingly we separate the velocity field 
into its mean and fluctuation, U = Vx + u. Below we use 
the wall units: y+ = y7fe T /L and V + (?/ + ) = V/ ^Jp'L. 

In the polymer laden flow, the mean velocity pro- 
file as a function of the distance from the wall (in wall 
units) exhibits an increase with respect to the Newtonian 
flow (see Fig. [1]). This is the phenomenon of drag reduc- 
tion. Note that the relative smallness of 7fe means here 
that the effect is not large, and by y + = 80 the velocity 
profile of the polymer laden flow is already parallel to 
that of the Newtonian flow. So any comparison with the 
theory of drag reduction should be limited to the rather 
narrow window of 20 < y + < 80. One should note the 
very different situations here and when 7& and Ve are 
very large. In the latter case, the mean velocity profile 
attains the maximum drag reduction asymptote (MDR) 
and never becomes parallel again to the von Karman log- 
law of the Newtonian flows. 

In Fig. O we show the momentum fluxes. There is 
a significant reduction in the Reynolds stress W(y) = 
— (u x u y ) for the rodlike polymer laden flow as compared 
to the Newtonian flow. The Reynolds stress is the mean 
mechanical momentum flux from the fluid to the wall. 
It had been explained before that the reduction in the 
momentum flux is at the heart of the mechanism for drag 
reduction [7j. Let us note that Figs. [1] and [2] are quite 
close to what one observes in turbulent channel flows with 
flexible polymers: an increase of the mean velocity profile 
and a marked decrease of the mean momentum fluxes. 



The theory of turbulent drag reduction by rodlike poly- 
mers [6] is based in part on the exact momentum balance 
equation: 

(<T xy )+is S + W = P '(L-y) (4) 

where S(y) = dV(y)/dy is the mean shear. A central 
ingredient of the theory is the statement that the polymer 
contribution to this equation, i.e. {(J xy ), can be evaluated 
as: 

(a xy ) « c 1 iy p R yy (y)S(y) (5) 

with some constant c\ and Rij = (TZij ) . Similarly, in the 
energy balance equation: 

vo(sijSij) + (crijSij) w SW (6) 

the polymer contribution to the dissipation, denoted here 
as e p = (cr^Sij), can be evaluated as 

e*>(y)*c 2 v p R yy (y)^- (7) 

with some constant C2, where K(y) = (\u\ 2 )/2 is the ki- 
netic energy of the fluctuating velocity. Here Sij is the 
fluctuating part of the velocity gradient tensor S a b, de- 
fined by: 

<S a6 (r, t) = S(y)S ax S by + s ab (r, t), (s ab (r, t)) = 0. (8) 

It was shown that for large 7^e and Ve these equations 
predict the establishment of a new velocity profile, again 
in the form a power law, but with a considerably larger 
slope compared to the Newtonian slope. This asymp- 
totic log-law is known as the Maximum Drag Reduction 
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FIG. 3: The profiles of the averages of the components of the 
conformation tensor R xx (solid line), R yy (circles), R zz (tri- 
angles), and R xy (squares). 



Asymptote (MDR). Moreover, the existence of a new log- 
law, with V + linear in log?/ + , for the (drag reduced) 
mean energy profile is directly related to R yy (y) increas- 
ing linearly with y. Physically, the theory states that the 
effects of the polymer can be treated as an ^/-dependent 
effective viscosity which increases linearly with y. It is 
thus of immediate interest to test these predictions also in 
the present case of relatively low 7& and Ve. To this end, 
we show in Fig. [3] the averages of various components of 
the conformation tensors: R xx , R yy , and R xy obtained 
in the simulation. We see clearly that R yy increases lin- 
early with y up to y^ ~ 80, which is the relevant range 
where drag reduction takes place in this simulation. In 
Fig. [U we present the direct comparison of (cr xy ) with 
v p R yy (y)S(y). The good agreement between the object 
and its evaluation is shown to exceed the region of linear- 
ity in Thus two central predictions of the theory are 
well supported by the direct numerical simulation even 
at the modest value of TZe that is available here. 

Needless to say, not every prediction of a theory that is 
developed as an asymptotic theory in the limit of Jte —> oo 
and Deborah number Ve — > oo [6] can be expected to hold 
verbatim, and some modification might be required. In 
the asymptotic theory one argues that R xx ~ 1 Rxyi 
R xy ^> R yy and R yy ~ R 2 y . We see from Fig, [3] that 
this is not the case here. To understand these results, 
note first that the Deborah number for flows with rodlike 
polymers is defined in the literature dLsVe^S/js, where 
7^ is the Brownian rotational frequency. The simulations 
were done using jb — 0, formally at infinite Ve for a lam- 
inar shear flow. In the present case of a turbulent channel 
flow, the effect of turbulence is to induce rotations of the 
polymers, giving rise to an effective relaxation frequency 
Iturb which depends on the turbulent intensity. Thus the 
corresponding effective Deborah number Ve = S/^turb is 
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FIG. 4: A comparison of R yy (y) (dashed line) with 
{(J xy ) /civ p S(y) (circles). 

also finite. At finite Tte and finite Ve, the relative sizes 
of the averages of the various components of the confor- 
mation tensors need to be reevaluated. In the following, 
we shall show how the theory can be employed for the 
case of finite 7^e and Ve to explain the observed numerical 
results Fig. [U 

We shall start from the equations of motions of the 
conformation tensor 7£ a &. Averaging Eq. ([2]) over the 
turbulent fluctuations, 

(9) 

where E = (TZ a bS a b)- To derive this equation we first 
employed the closure assumption 

{HabcdScd) ~ {TtabUcdScd) « R a b{RxyS + £) • (10) 

The second simplification is the assumption that after 
removing the mean shear, the remaining velocity fluctu- 
ations are not too far from isotropic, and in the log-layer 
can also be taken as homogeneous. This implies that cor- 
relation functions of with lZ a b are isotropic in space: 

(Sxc'R'Cx) ~ (Syc'R'cy) ~ zc^cz) "~ (H) 

J] 

(SacKcb) ~ $ab— (12) 

Within the same assumptions we can also offer an ap- 
proximate evaluation of E w A^/K/y 2 since all the 
velocity fluctuations are close to isotropic. Here A is 
a constant of the order of unity. Finally, we evaluate 
{UkdlZab/ dxk) ~ 0. This is seen by integrating by parts 
and using the the incompressibility constraint on. The 
derivative of the average is negligible for fluctuations that 
are not too far from homogeneous. With all these we ob- 
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FIG. 5: Theoretical prediction of the averages of the compo- 
nents of the conformation tensor as a function of the distance 
from the wall. 



FIG. 6: Measured (solid line) and predicted (dotted line) dis- 
sipation contributed by the rodlike polymers to the energy 
balance equation. The measured quantity e p has two peaks, 
whereas R yy K/y 2 has only one peak. 
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These equations are identical to those obeyed by the con- 
formation tensor in a steady laminar shear flow and this 
tells us that S/H can be taken as the effective y- dependent 
Deborah number in the simulation. Using Eqs. (fT3]) - (fT5]h 
we proceed to compute the y profiles of R^. To do so, 
we solve R^ in terms of Then we consider the 

polymers to be small perturbations and use the momen- 
tum and energy balance equations for turbulent Newto- 
nian channel flow to get S(y) and K(y). Written in wall 
units, the momentum and energy balance equations for 
turbulent Newtonian channel flow read 



S+ + W + = 1 - -s— 



Re T 



3/2 



W+S+ . 



(16) 
(17) 



where is the Von Karman constant and 5 is the thick- 
ness of the viscous layer, and S « 6 was found in 0. 

In Fig. [5l we show R xx , R xy and R yy obtained by solv- 
ing $131), CHJ) and (US]) with A = 0.7. We find general 
agreement with the results shown in Fig. [3l explaining 
why in this case the relative sizes of the averages of the 
components of the conformation tensor differ from the 
predictions of the asymptotic theory. Our simple mod- 
eling can also explain the results shown in Fig. 2) We 
obtain from Eq. ((Til) 



SRyy 2R X y(S R X y + X) 



(18) 



Using this result and Eqs. (J3j) and ([TO)), we get 



: 6VpR X y(R X yS - 



QlSpRyy S 



(19) 



Thus, even if the asymptotic theory cannot be applied 
directly in the present DNS, the basic prediction a xy « 
c\v p RyyS still holds, consistent with the above men- 
tioned approximations, i.e. neglecting anisotropic con- 
tributions like (R xa Say)' 

In order to state that R yy is an effective viscosity it 
should also play the role of additional viscosity in the 
energy balance equation. To test the validity of this 
we compare the energy dissipation due to polymers e p 
against the theoretical prediction Eq. 0. This is done in 
Fig.[6l On the one hand, the results shown in Fig. [6] sup- 
ports our general conclusion that e p ~ R yy K/y 2 within 
a prefactor that we cannot estimate from the theory. On 
the other hand the agreement is still not perfect; the mea- 
sured e p shows two maxima as a function of y + . The first 
peak is related to the maximum of R yy K/y 2 while the 
second peak corresponds to the maximum of (R xx s xx ). 
One of the central statements of the asymptotic theory 
@ is that this term drops, by exact cancellation with an- 
other term, when the MDR is approached in Ve — > oo. 
This is not occurring yet in our simulations with low 
7^e and Ve. 

From our numerical simulations and the theoretical 
analysis we can also state that R yy is linearly growing 
up to y + ~ 80. To see this, estimate the effect of linear 
viscosity profile using the following equations @ : 



[1 + a(y - 5)}S+ + W+ = 1 , 



y 



— = S+ . 



(20) 
(21) 



In asymptotic conditions F = 1, but here the factor F 
takes into account that the effective slope of the linear 



5 




FIG. 7: Comparison of the profiles V + {y + ) and W + (y + ) as 
computed from the balance equations ([20]) and ([2T|) with the 
results of numerical simulations. The agreement is excellent. 



viscosity profile is somehow smaller for the energy bal- 
ance equation than for the momentum equation. Ac- 
cording to our previous discussion on Fig. ([6]), we can 
estimate F « 2. In (|2T|) the term A(a) was determined 
theoretically Q: 



A(a) = 



1 — aS 



In our case the value of a is given by the relation 

Vp (IRyy 

v dy 



(22) 



(23) 



where the slope dR yy /dy is estimated from the numerical 
simulation. Solving these simple coupled equations we 
present in Fig. ([7j) the quantities V + (?/ + ) and W + (?/ + ) 
respectively for the Newtonian flow and for the polymer 
laden flow. Both figures agree well with the DNS. 



CONCLUSIONS 

The general philosophy behind our approach to drag 
reduction by additives is to consider the balance equa- 
tions for mechanical momentum and turbulent energy, 



and to analyze the predictions of these equations for the 
profiles of the relevant quantities, in particular the mean 
velocity at distance y from the wall Q . The theory can be 
simplified in the asymptotic regime when TZe and Ve are 
very large; there one finds universal profiles, in particu- 
lar for the the mean velocity profile which becomes the 
universal MDR [9[. For comparison with numerical sim- 
ulations, where the drag reduction effect is rather limited 
due to small finite 7^eand Ve, one needs to analyze the 
balance equations with greater care, taking into account 
the non-asymptotic effects. This is what we have done 
in the present work, and our results are shown to agree 
well with the data obtained from direct numerical simu- 
lations. 
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